Generic flow proflles induced by a beating cilium 
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We describe a multipole expansion for the low Reynolds number fluid flows generated by a localized 
source embedded in a plane with a no-slip boundary condition. It contains 3 independent terms 
that fall quadrat ically with the distance and 6 terms that fall with the third power. Within this 
framework we discuss the flows induced by a beating cilium described in different ways: a small 
particle circling on an elliptical trajectory, a thin rod and a general ciliary beating pattern. We 
identify the flow modes present based on the symmetry properties of the ciliary beat. 
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I. INTRODUCTION 



Cilia are thin cellular protrusions that beat in an asym- 
metricperiodic fashion in order to propel the surrounding 
fluid [ll-Q. They are involved in swimming and feeding 
of a number of protozoa and also have many crucial func- 
tions in vertebrates. These include the mucous clearance 
in respiratory pathways, transport of an egg cell in Fal- 
lopian tubes, left-right symmetry breaking in embryonic 
development 0-0, and, recently discovered, otolith for- 
mation in hearing organs Q . Cilia have inspired designs 
for microfluidic pumps and mixers using similar beating 
patterns [Ml. 



The beating pattern of a single cilium can be very com- 
plex and generally consists of a working stroke during 
which the stretched cilium moves the fluid in the direc- 
tion of pumping and a recovery stroke during which the 
cilium folds and returns to the origin by sweeping along 
the surface, thereby minimizing the drag as well as back- 
ward flow. Especially when working in large ensembles, 
the beating patterns found in nature are remarkably close 
to the theoretically calculated optimum [15]. Despite this 
complexity, the flows can be described with a small num- 
ber of generic terms in the far field, i.e., when we observe 
the flow at a distance r which is sufficiently larger than 
the cilium length L. In the leading order the flow velocity 
falls as 0(r-2) [l^. 



II. MULTIPOLE EXPANSION 

A. In unbounded fluid 

In the low Reynolds number regime, the motion of an 
incompressible fluid is described by the Stokes equation 
and the incompressibility condition: 



?7Av — Vp 
V-v = 0. 



(1) 
(2) 



To derive the general form of the flow induced by a local- 
ized distribution of forces and/or sources we follow the 
approach of Lamb [l^ , also used by Happel and Brenner 
flSj . Note that this is just one of many complete sets of 
far-field solutions (see Ref. [19] for a derivation in Carte- 
sian coordinates). Although the basis set we will de- 
rive can be obtained more directly as spatial derivatives 
of known solutions (Stokeslet and source), the following 
construction has several advantages. First, we know from 
the beginning on that we are dealing with a complete set 
of linearly independent solutions. Second, these solutions 
will appear classified by their angular symmetry, which 
are useful for the description of flows induced by cilia 
which also share some of these symmetries. And finally, 
the construction allows us to distinguish between differ- 
ent solutions (e.g., with or without a pressure gradient). 

Multiplying Eq. ([T]) with V- immediately leads to a 
Laplace equation for the pressure 



Aj9 = . 



(3) 



In this paper we go beyond the quadratic order and 
provide a complete set of functions that describe the flow 
in the presence of a no-slip boundary. We then specifi- 
cally discuss the 0{r~^) terms that are present with dif- 
ferent beating patterns, depending on their symmetry 
properties. 
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The solution of the homogeneous equation in the absence 
of pressure gradients (Av = 0) can be written as v = 
V X (rx) + V^, where x and <l> are two solutions of the 
Laplace equation. In total we need three sets of harmonic 
functions (p, x and ^) to construct the general solution 
of the Stokes equation in unbounded space. We now 
introduce a spherical coordinate system such that 



(4) 
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If there is no external flow, i.e, the hm^^oo v(x) = 0, the 
multipole expansion of the velocity fleld reads 



oo 

E 

1=0 



vi with ~ r" 



(5) 



with 



l-l 



27^(1 + 0(1 + 20 



2 + / 



7/(1 + 0(1 + 20 



and 



^-1 - i 



^'i-i(cose) 



(6) 

(7) 
(8) 



Pl+^ = 7?(1 + 0(1 + 20+2 E C'ime'"**P,™i(cos^) . 
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The coefficients Aim are deflned for / > 1 and |7ri| < 
Bim are deflned for / > 1 and |777.| < / — 1 and Cim for / > 
and 1 772 1 < / + 1. Using an elementary vector identity the 
flrst term in (|6j) can alternatively be written as (VxO x r. 
Inserting these terms into Eq. (j6j) gives the expression for 
the velocity 
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In the unbounded space, the leading terms have the or- 
der / = and represent the 3 components of a Stokeslet. 
Their magnitude is given by the coeflicients Co,-i, Co,o 
and Co,i. In the next order, / = 1, representing terms 
that decay as 1/r^, we have 3 solutions for the x com- 
ponent, 1 for ^ and 5 for p, 9 in total. In general, the 
number of terms of order / is 



iVz = 3 + 6/ . 



(11) 



Note that all these solutions can be constructed from 
derivatives of 4 fundamental solutions. These include 3 
Stokeslets (/ = 0, 77i = 0, ±1) and a source / = 1, 77i = 0. 
We denote these solutions as Q^^ Qy, Qz and S. The 
remaining solutions can be constructed from derivatives 
like dj^Qx, dyQ^ etc. 



B. Bounded by a no-slip plane 

As the cilium (or any other source of fluid pumping) is 
embedded in a plane, we have to flnd the solutions that 
fulflll the boundary condition 



v(r,6> = 7r/2,(^) =0 



(12) 



We flrst note that this condition has to be satisfled for 
all r and (j) values, which means that we can collect terms 
with the same indices / and m and then form all indepen- 
dent linear combinations that fulflll the boundary condi- 
tion. 

When expressing vim in the plane (at ^ = 7r/2) we have 
to distinguish between even and odd values of / + 772. If 
/ + 771 is even, so is the associated Legendre polynomial 
P[^{x), implying -£^Pl^{x)\^^^ = 0. For an odd / + rri, 
PJ^ix) is an odd function so that P^^(O) = 0. Of course, 
we have to take into account that the terms with an even 
/ + 771 actually contain $ and p terms with an odd l^m, 
because they contain spherical harmonics of the order 
/ — 1 and / + 1, respectively. In total the condition that 
vzm(^, 7r/2, ^) = for an even / + m implies 



imAimPr{0)-BimP!!I',{0) 



CimPl^,{0) = 



(13) 

This condition is generally fulfllled by 2 independent so- 
lutions. For example, we can choose the coefficients Aim 
and Cim freely, but then have to determine Bim from 
(p!3|) . However, when m = ±1, there is a single solution, 
because there is no ^ {B) term. When / = 772 = there 
is neither a ^ nor a x term, so there is no solution of this 
order. 

For an odd / + 777, the boundary condition leads to two 
equations 

Ai^PnO) + imBi^Pi"^,{0) + im^CimPl^li = 

-ZP,™P,™i(0) + ^ ^ ^ Ci^PlliiQ) = 

(14) 

These equations generally have a single solution for each 
pair 777. 

It is now time to count the total number of independent 
solutions for each /. An overview is given in Table HI For 
any / > 1, we have 

Ni = 31 (15) 

solutions, down from 3 + 6/ in the unbounded space. 

In the following we will have a closer look at the terms 
that decay quadratically with the distance (/ = 1) and 
those that decay with the third power (/ = 2). 



C. 1 = 1: Modes of the order r 

In this order we have Ni = 3 modes: one of them 
(777 = 0) represents a fluid source and two {m = ±1) the 



3 



TABLE I: Number of linearly independent solutions for each 
pair of numbers / and m. 



m ±1 


±2 


±3 ±4 . . 


. Total 


/ = 









1 1 1 






3 


2 2 1 


1 




6 


3 1 2 


1 


1 


9 


4 2 1 


2 


1 1 


12 


I 






31 




flow of two horizontal Stokeslets in the proximity of the 
boundary. 

The source is described by the second and the fourth 
term in Eq. ([6]) alone. The solution for the pressure is 
Pio 0^ _^ 3cos^ 6>-i ^ ^^^^ second term <l>oo cx 1/r. 

The resulting velocity is 



S'(x) = — cos 



y 



(16) 



To distinguish it from a source in bulk we call this term 
^'surface sourcelet". Note that because of volume conser- 
vation beating cilia naturally do not act as a source. But 
the source is one of the fundamental singularities and we 
will later show that its derivatives can be used to express 
higher order terms that are present. Also, this formalism 
could be used in a more general context - the source could 
describe fluid injection through a pore in the surface, or, 
temporarily, the flow around an expanding bubble. 

The solutions with m = ±1 represent two horizon- 
tal Stokeslets. In this order there is no ^ term and the 
X term has to vanish {Aim = 0) in order to fulfill the 
boundary condition (p!3|) . The velocities are determined 
by the last term in Eq. (|6]) and stem from pressure terms 
P2,±i ^e^*^"^ cos^sin^. We can construct the first 
from vi^i + vi^_i 



COS U sm U cos (per = 



(17) 





and the second from vi i — vi _i 



Dy = — cos i 



yz 

sin 6> singer = —r 



(18) 



FIG. 1: Fundamental modes of the order . For each mode 
the angular dependence of the velocity field is shown along 
a semi-transparent hemisphere. The colour of the arrows de- 
notes inbound and outbound flow. The three depicted modes 
are a "surface sourcelet" S and two "surface Stokeslets" Dx 
and Dy. 



We will call these terms ^'surface Stokeslets^\ With 3 
linearly independent solutions we know that we have a 
complete set for the order v cx r~^. 

The surface sourcelet 5* and the two surface Stokeslets 
Dx and Dy are shown in Figure [TJ Note that these 3 



modes represent a fundamental set. All higher order 
terms can be constructed by deriving them by x and y. 
We can see this from the fact that we have counted 3/ so- 
lutions of order /. To get from the fundamental solutions 
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to those of order / we have to derive them / — 1 times 
by the coordinates, which gives us / hnearly independent 
derivatives (unhke in unbounded fluid they ah are inde- 
pendent). So we have shown that in total 3/ derivatives 
of the fundamental solutions represent a complete set of 
solutions to that order. We will nevertheless derive the 
solutions from the spherical harmonics directly in order 
to collect modes according to their angular symmetry - 
this system is very convenient to discuss the ciliary flows. 



D. 1 = 2: Modes of the order r~'^ 

In this order we have a total of 6 independent terms, 
two with m = 0, two with m = ±1 and two with m = ±2. 
The first solution can be obtained from X2,o ^ 

1 3cos^g-l . 



Qi = ^ sin ^ COS ^ 



z 



l-y 

X 

\ 



(19) 



This is a rotlet around the z axis. 

The second solution is obtained from i>i_o and pa^o- It 
reads 

1 



Q2 = ^ ((4 cos^ 6* - 2 cos e)er 



z 
I7 



/ (5^2 
(5^2 



2r2)y 



(20) 



and corresponds to the flow of a vertical Stokeslet near 
the boundary plane [20]. The fluid is moved outwards 
along the z axis and inwards along the x — y plane. Note 
that this is the only axisymmetric term - it is there- 
fore clear that this is the leading term describing acous- 
tic streaming caused by ultrasonic oscillations of a small 
bubble on a planar surface [2l|. 

The modes with mil that fulfill the boundary con- 
dition (p!4|) are 



Qz = (4 cos^ d sin Qer - cos^ 




icos^ l9e^) e^^ 



(54 = ^-3 (4 cos^ d sin — cos^ deQ — i cos^ 



(21) 



(22) 



They represent the fields of a horizontal rotlet (around 
the y and x axes) in the presence of a boundary. As we 
will see later, they can also be viewed as source-sink pairs 
(source-doublets) on the surface. 



Finally, the modes with m = ±2, which have to fulfill 
the boundary condition ([13]) are 



5R— sin cos (4 sin Ocr — cos i 

( (5(x2-^2)_^2)^\ 

(5(x^ — y'^) + P)y 
\ 5{x^-y^)z J 

sin 6 cos 6 (4 sin Ocr — cos i 
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(23) 
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( (10x2 - r^)y " 
(10^2 — r'^)x 
y lOxyz 



(24) 



These modes represent surface stresslets. Q5 is generated 
if a pair of forces pulls the fluid apart along the x axis 
and another pair together along the y axis. Qq is similar, 
but rotated by 45° about the z axis. 

Like in unbounded space fioj all these modes can be 
represented by derivatives of the ~ r"^ modes in the 
following way: 



Qi 
Q2 



Qg 
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dy 


dx 
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dD, 


dDy 


dy 


dx 



(25) 
(26) 
(27) 
(28) 
(29) 
(30) 



These relations indicate that the modes Qi, Q2, Qb and 
Qq can be seen each as the flow of four point forces at 
an infinitesimal distance, as well as infinitesimally close 
to the surface. Qa and Q/^ represent a dipole consisting 
of a source and a sink. This representation is illustrated 
in Fig. [31 



E. Higher order modes 

In each higher order /, we have A^^ = 3/ modes that 
fall off with the power ^-(^+1). All of them can be con- 
structed from the derivatives 



dx^dy 



i-i- 



dx^dy^~^~ 



and 



dx^dy^~^~^ 

(31) 



Unlike in unbounded fluid all these derivatives are lin- 
early independent. 



Qi Q2 




Q5 Qe 



FIG. 2: Six modes of the order r ^: The surface rotlet Qi, the vertical Stokeslet Q2, the surface source doublets Q3 and Q4 
and surface stresslets Q5 and Qe- 

F. Normal form of the flow proflle To next order, the flow now has the shape 

Ciliary flows are characterized by the fact that they 
do not contain a source, but they do generate directional 
flow, so at least one of the terms and Dy is present. 
By rotating the coordinate system, we can always bring 
the time-averaged flow to the form v(x) = ADx-\-0{r~^). 



6 

v(x) = AD,{^) + ^S,A(x) + 0{r-^) . (32) 

i=l 

If we transform the coordinate system to = :>c — dxCx — 



6 



'J 



Q4 



^^^^^ 



05 



Qs 



y 



FIG. 3: Six modes of the order r ^ represented with funda- 
mental solutions. 



dyCy, we have 



V'(X') = V(X' - dxCx - dyCy) 

= ^(D.(x')-d.^(x')-d,^(x') 

6 

+ ^BiA(x') + 0(r-^) . 

i=l 

As we know from Eqns. (|25ti30p the derivatives are 



dx 



dy 



(33) 



(34) 



and by choosing dx = — and dy = — we can elim- 
inate the modes and Qq. Of course, this transforma- 
tion is only possible if ^ 7^ 0. Thus we have brought the 
stationary flow to the normal form 

4 

v'(x') = ^D,(x') + ^ i3pi(x') + 0(r-4) . (35) 

i=l 

To summarize, with the proper choice of the coordi- 
nate system the static flow generated by a cilium can be 




X 



^^^^^^^^ 



X 



^^^^^^^^^ 



FIG. 4: a) The volume flow rate Q is defined as the fluid flux 
through a half-plane perpendicular to the surface and to the 
direction of pumping, b) Above a densely ciliated surface the 
flow becomes homogeneous with the average velocity Vc- 



described up to the order with only 5 components: 
the surface Stokeslet Dx, the surface rotlet Qi, the verti- 
cal surface Stokeslet Q2, and the surface source-doublets 
Qs and Q4. Because we can only transform the coordi- 
nate system once, for the static flow, we cannot reduce 
the oscillatory flows in the same way. 



G. Volume flow rate 

The performance of a cilium is characterized by the 
volume flow rate Q, defined as the average flux through a 
half-plane perpendicular to the direction of pumping 22, 
Fig. |4^. For example, for the mode corresponding to the 
pumping in x-direction with amplitude A, v(x, z) = 
ADx{x, y, z) has the volume flow rate 



.00 .00 2 
Q=l dy I dzvx{x,y,z) = -A . 
J-00 Jo 



(36) 
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Due to volume conservation the flow rate is, of course, 
independent of the x-position of the cross-section chosen. 
Similarly, the volume flow rate of the mode Dy would be 
deflned as the flux through a half-plane with a constant 

y- 

For the source, the flow rate is the actual influx. For 
a flow with v(x, z) = AS{x^ z)^ the flow rate is 



where R = (X, F, Z) is the position of the particle, Z is 
its height above the plane, R = (X, F, —Z) is the position 
of its mirror image and 



(37) 



The volume flow rate for higher order modes (/ > 2) is 
always zero - the easiest way to see this is from the fact 
that they can be represented as spatial derivatives of the 
fundamental modes. 



H. Velocity above ciliated layer 

Instead of a single cilium, we often deal with a densely 
ciliated surface (Fig. HJd). Let p denote the surface den- 
sity of those cilia. Then the velocity above an infinite 
array is 



dX 



27T 

dYpADx{x-X, y-Y, z) — pA = irpQ 



(38) 

which is independent of z. The second relationship ex- 
presses the velocity above a ciliated layer with the volume 
generated by each cilium. In this regime, one can simplify 
the description of cilia by replacing them with a surface 
slip term with velocity Vc [23] . 



III. CILIARY MODELS 

The ciliary beating pattern is quite complex, but there 
are several simplified models that have the correct sym- 
metry properties and capture the far-field fiow patterns. 
The simplest of all is a model that replaces the cilium 
with a single small particle (radius a), moving periodi- 
cally on a closed path, and was used in several models 
for ciliary synchronization ■24 - 27J . 



87rr/ \ |x| Ixf 



and 



Gf^(x) = (l-25^,)^GL(x) 



(40) 
(41) 

(42) 



are the fields of a Stokeslet, source doublet and a 
Stokeslet-doublet, respectively. A Taylor expansion of 
Blake's expression for small particle located at (0, 0, Z) 
up to the order reads 

v(x) = ^(l2Z(F^I)^(x) + F^I)^(x)) 

+ 6Z2(-F,Q3(x)-F^Q4(x) + F,Q2(x))) . (43) 

As noted by Blake and Chwang [20] the horizontal forces 
Fx and Fy generate terms of the order r~^, while the 
'vertical force Fz only generates terms ~ . If the par- 
ticle is not located exactly on the z axis, but in any point 
(X, F, Z), we can use the derivatives (|25H3Q|) to obtain an 
expression that is again exact up to the order r~^: 

v(x) = ^{l2Z{FxDx{-s) ^ FyDy{-K)) 

+ 6Z2(-F,Q3(x) - F^Q4(x) + F,Q2(x)) 

+ 6zxF,(g5(x)-g2(x)) 
+ 6zrF,(g6(x)-Qi(x)) 

+ 6ZXF^(Q6(x) + Qi(x)) 

+ 6ZrF^(-Q5(x)-Q2(x))) . (44) 

This equation provides a basis for discussing all models 
that replace the cilium with a point particle. 



1. Oscillating particle 



A. Cilium as a small sphere 

We start our discussion by expanding the fiow field of 
a small particle subject to force F in terms of D- and Q- 
modes. An exact solution for the fiow has been derived 
by Blake [28] and is now sometimes referred to as Blake's 
tensor: 

v«(x) = ^[Gf^(x-R)-Gf^(x-R)+2Z2G^^(x-R) 

-2ZGf^(x-R)]F^ (39) 



We will start our discussion with a small particle, 
whose position is oscillating periodically parallel to the 
X-axis (Fig. [5^). Of course, this reciprocal motion will 
not generate any net (time-averaged) fiow [29]. However, 
this type of motion has been studied as a model system 
for ciliary synchronization [g^, |30| and therefore its far- 
field is of interest. It should also be noted that even if a 
single cilium performs reciprocal motion, an array of cilia 
with metachronal coordination can still generate directed 
fiow [31, 32]. 

Since the motion is symmetric with respect to the x — z 
plane, the only solutions that are theoretically possible 
are D^, Q2, Qs and Qs- 
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d) e) f) 

FIG. 5: Models describing the ciliary beat as the motion of a small sphere along a closed trajectory: a) Linear, back- and- forth 
movement, which does not generate a directed flow, b) Planar elliptical trajectory, c) Tilted elliptical trajectory, d) General 
trajectory shape, e) General shape at a constant height, f) A model describing the cilium as a stiff slender rod. 



We parameterize the particle motion as 



R(t) 



(b 



V 



sincjt 

D 



(45) 



so that the force acting on it is F = QirrjaB cos {ujt)ex. 
The resulting velocity field is 

v(x,t) = 5^^^^ [sm{ujt) (12I)^(x) - 6DQs{^)) 



and the force on the particle is F = 67rr]a{B cos{ujt)ex — 
C sm{ujt)e^ . 

Based on the symmetry of the trajectory, which is in- 
variant when we transform x ^ —x and t ^ —t simul- 
taneously, we conclude that all terms that are even in t 
have to be even in x. That reduces the possible terms to 
Dx and Qs- Conversely, the terms that are odd in time 
have to be odd in x, which restricts them to Q2 and Qs- 
The velocity field up to the frequency u (we omit higher 
harmonics, which exist up to the frequency 3uj) reads 



sin(2c^t)35(Q2(x) - Q,{^))] . (46) ^(^^ t) = — [6BCDx{^) - 6BCDQs{^) 



2. Planar beat 

The next level of model is a planar ciliary beat. We 
model it with a small particle moving along an elliptical 
trajectory in the x — z plane (Fig. [Sb). Individually, 
planar cilia are rather inefficient [33], but nearly planar 
beating patterns can be found in respiratory epithelia, 
as well as in some microorganisms, such as Opalina [2|. 
Several models for ciliary synchronization study planar 
beating patterns [33-[36j. Because of the symmetry with 
respect to the x — z plane, the solution will still consist 
of the following terms alone: Dx^ Q2, Qs and Q5. The 
particle trajectory is now parameterized as 



{l2BDDx(:x) - 6BD'^Qs{x.)) cos(cjt) 

+ Ib^C{Q,{^) - g2(x)) sm{ujt)] . (48) 



3. Tilted ellipse 

More properties of a ciliary beat are captured by a 
model that describes it as a sphere, circling on a tilted 
elliptical trajectory as shown in Fig. [5fc. It is de- 
scribed as 



( B 



sincjt 




^ Bsiiiujt ^ 

Acosujt 
\C cosuot^ D J 



(49) 



(47) 



\C cosujt ^ D J 



The {x ^ —x^ t ^ —t) symmetry is still present, so that 
the stationary current can only consist of components 
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Dx, Qi, Qs and Qq. It is 



v(x) 



3(jua 



[6BCDx {x)-6BCDQ3 (x) -6 A5L>Qi (x) 



(50) 



The major difference between the tilted ellipse and the 
vertical ellipse is that the tilted one contains a Qi com- 
ponent, or a surface rotlet, as an obvious consequence of 
the rotation around the z axis. In the case of a flat ellipse 
{C = 0), this is the only stationary component. 



4' General trajectory 

Let us finally discuss the stationary velocity field of a 
small particle moving along an arbitrary periodic trajec- 
tory (Fig. [SJi). According to Eq. dH]) the velocity is 



v(x) = Y Dx{x.) ^ ZdX^Dy{x.) j ZdY 

- iQ3(x) j Z^dX - iQ4(x) ^ Z^dY 

9M_Q^ fzxdx^9M_QMfzYdx 



+ 



Q6(x) +Qi(x) 



(51) 



where T denotes the period. The first two terms (through 
a coordinate rotation they can be reduced to just Dx) 
generally dominate and they show that the volume flow 
rate is proportional to the projection of the trajectory 
onto the vertical plane in the direction of pumping. An 
interesting special case arises when the particle does not 
move vertically, Z = const (Fig. [5^). Then all integrals 
except those with YdX and XdY represent integration 
of a total differential along a closed path and therefore 
vanish. Equation (|5T]) simplifies to 



v(x) 



9aZS 



Qi(x) 



(52) 



where 5* is the area of the trajectory, S = § XdY = 
— § YdX. Regardless of the exact shape of the trajectory, 
a particle moving in the horizontal plane only produces 
a surface rotlet in the order . 



B. Stiff slender rod 

The next level of complexity is to describe the cilium as 
a thin stiff rod. Although the hydrodynamic description 
in this model is more accurate than replacing the cilium 
with a small particle, the motion is even more restricted. 
At any time it is determined by two angles. There are 
different ways to calculate the forces on the cilium. The 
simplest is the local drag or resistive force theory (RFT) 



which assumes two constant friction coefficients for mo- 
tion parallel and perpendicular to the orientation of the 
rod [33] • It has the advantage that it allows analyti- 
cal calculation of the force distribution, while still giving 
good results. Alternatively, one could use the slender 
body theory (SBT), which includes hydrodynamic inter- 
actions for distances longer than a chosen cut-off |38| . 
Even better accuracy can be achieved by describing the 
cihum as a chain of spheres 0, ^ , a chain of regularized 
Stokeslets ^] or with surface boundary elements [40]. 

Within the resistive force theory the local drag is de- 
scribed by a tangential drag coefficient Ct and a normal 
coefficient Cat - in the case of a rod that is anchored with 
one end only the latter is relevant. From the force den- 
sity, we can calculate the far-field fluid flow using Blake's 
tensor ([39|) or its far field approximation (|44|) . 

A simple version of the slender rod model is obtained 
if the rod moves along the mantle of a tilted cone, as 
proposed by Smith and coworkers [22]. The motion of a 
point on the cilium is then described as (Fig. [5]F) 



R 



cos sin ip cos ujt — sin 6 cos ip 
y sin 9 sin ^jj cos uot + cos 9 cos 



3L 



V 



/ 

sin 9 cos ip 


(53) 

with s denoting the position parameter running from 
to L. We have placed the cone such that its center of 
weight, not the tip, lies on the z axis. We will soon 
see that this way we eliminated the terms Qs and Qq 
and brought the flows to the normal form discussed in 
Section III F[ In RFT the force density on the cilium is 
simply the transverse drag coefficient, multiplied by the 
local velocity 

/ sin ip cos ujt \ 



f = Cn (^s 



(54) 



- COS ^ sin sin cjt 
y — sin 9 s'milj sin ut J 

We obtain the far-field flow by inserting the force densi- 
ties (|5^ and positions (|53]) into Eq. (|^^ and integrating 
over the length. The resulting average flow is 

v(x) = sin^ ip (sin 9Dx — -rL cos^ 9 cos tljQi 



— -Lsin6>cos6>cos?/^(53j . (55) 

The first term is the well known surface Stokeslet term 
Dx, with an amplitude proportional to sin 9 sin^ ip [22], or 
the projection of the tip trajectory onto the x — z plane. 
It is maximal for = arctanV^ and 9 = it /2 — ip. The 
second term, Qi, is a surface rotlet and is maximal for 
^ = (no tilt) and ip = arctan ^/2. The third term, Qa, a 
surface source-doublet, appears almost inevitably along 
with the Dx term. 

An approximate value of the friction coefficient Cn was 
given by Gueron and Liron 

= l + 2H2q/a) ^^^^ 
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where a is the radius of the rod and q some character- 
istic length scale of the order \fLa. A comparison with 
the Rotne-Prager approximation gives a more precise es- 
timate of Cn = 1.22717] for the ratio L/a = 20 [42]. 



C. General 

In the following we will discuss the flows generated 
by an arbitrary periodic beating pattern, based on its 
symmetry properties. Formally the problem could be 
solved numerically within the resistive force theory (43| . 
as described in the previous section, although the motion 
would no longer be purely transversal, so the tangential 
friction coeflicient would need to be included, too. Of 
course, any of the more accurate methods mentioned in 
the previous section can be applied. 

A common symmetry of the beating pattern includes 
mirroring in the x direction and simultaneous time rever- 
sal {x —X, t — t), Fig.[6t. The patterns discussed in 
Sections IIII A 2\ IlII A 31 and IIII Bl all contain this symme- 
try. But note that many ciliary strokes found in nature 
do not. The typical recovery stroke during which the 
cilium bends and sweeps along the surface is not sym- 
metric. Interestingly, the theoretically optimal solutions 
sometimes break this symmetry spontaneously and some- 
times not, depending on the allowed radius of curvature 

The symmetry property of the cilium trajectory has to 
be reflected in the stationary flow component 



v(Tx) = 5Tv(x) 



with s 



-1 and 



T 




(57) 



(58) 



This condition is fulfllled by the modes D^, Qi, Qs and 
Qq. Note that this holds for the time- averaged flow (sta- 
tionary component). The oscillatory flow components 
can contain other modes. 

Another interesting symmetric case is the general pla- 
nar beating pattern (Fig. [6)3). Then Eq. (|57|) holds with 
5 = 1 and 



T = 




(59) 



It is fulfllled by modes D^^ Qs and Q5. An example 
of a planar beating pattern is a cilium producing a waving 
pattern, similar to a planar flagellum attached with one 
end to the surface. Opalina [2] produces this type of 
waves. 

A further special case contains planar patterns that are 
additionally symmetric with respect to the transforma- 
tion X -X, t ^ t ^ T/2 (Fig. [6t). This imposes the 



additional symmetry with 5 = 1 and 



T 




(60) 



which is only fulfllled by the models Q2 and Q^. An 
example is a planar flagellum that beats in a way that is 
symmetric in x direction, such that the waves propagate 
vertically. 

Let us flnally look at a rotationally symmetric beating 
pattern (Fig. [6]i). The symmetry transformation reads 
X Jx, t ^ t-^ At with 



^ cosip simp 



T = 



V 



- sin (f cos (f 
1 



(61) 



and the two modes that fulflll this symmetry are Qi and 
Q2. The flrst is the obvious surface rotlet. Interestingly, a 
rotating flagellum can also generate a vertical flow - this 
is in contrast with a single circling particle (|52]) . This 
type of flow is likely involved in the otolith formation in 
the developing ear of the zebraflsh embryo [13, EH . 



IV. CONCLUSIONS 

Although the ciliary beating patterns can be very com- 
plex and is not yet well understood from the mechanical 
perspective, their far-fleld properties up to the order 
can be described with only a few terms in a multipole ex- 
pansion. Of those, only the terms proportional to are 
fundamental - all higher order terms can be represented 
as their spatial derivatives. Depending on the symme- 
tries present in the stroke pattern, the number of terms 
describing the flow can often be further reduced. We 
have derived the magnitudes of different terms for sev- 
eral models describing the cilium as a small sphere, as 
well as for a rod-like cilium within the framework of the 
resistive force theory. Of course, the multipole expansion 
is not limited to any speciflc hydrodynamic approxima- 
tion, but more accurate descriptions require a numerical 
calculation of the amplitudes. 

Although explicitly applied to cilia, our approach is 
equally suited for any localized sources of flow in the 
proximity of a planar no-slip boundary. Examples in- 
clude an oscillating bubble [21j, a tumbling object [46], 
bacteria attached to the surface [i^ or swimming close 
to the surface [48] and many more. 

Let us flnally remark that our approach is not re- 
stricted to planar boundaries. Exact solutions exist for 
flows due to point forces in the vicinity of a sphere or in- 
side a spherical cavity [49] . Knowing the fundamental so- 
lutions on a sphere (surface source and surface Stokeslet), 
one can derive the higher order terms in a way similar 
to Eqns. (|25l[3Q|) . The derivation of the flows caused by 
a single cilium of a spherical swimmer (e.g., Volvox jsoj ) 
to the same order as shown here is then straightforward. 
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c) d) 

FIG. 6: Beating patterns according to their symmetry properties. In each figure the faint hne shows the cihum half a period 
later, a) A pattern symmetric upon mirroring in x direction and simultaneous time reversal, b) A planar beating pattern, c) 
A planar pattern that is symmetric with respect to x ^ — x, t ^ t + T/2. d) A pattern symmetric with respect to rotation 
around z axis and simultaneous time shift. This corresponds to a cilium rotating around the z axis. 
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